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ABSTRACT 

We study the formation of stellar haloes in three Milky Way-mass galaxies using 
cosmological smoothed particle hydrodynamics simulations, focusing on the subset 
of halo stars that form in situ , as opposed to those accreted from satellites. In situ 
stars in our simulations dominate the stellar halo out to 20 kpc and account for 30-40 
per cent of its total mass. We separate in situ halo stars into three straightforward, 
physically distinct categories according to their origin: stars scattered from the disc 
of the main galaxy (‘heated disc’), stars formed from gas smoothly accreted on to 
the halo (‘smooth’ gas) and stars formed in streams of gas stripped from infalling 
satellites (‘stripped’ gas). We find that most belong to the stripped gas category. 
Those originating in smooth gas outside the disc tend to form at the same time and 
place as the stripped-gas population, suggesting that their formation is associated with 
the same gas-rich accretion events. The scattered disc star contribution is negligible 
overall but significant in the Solar neighbourhood, where > 90 per cent of stars on 
eccentric orbits once belonged to the disc. However, the distinction between halo and 
thick disc in this region is highly ambiguous. The chemical and kinematic properties of 
the different components are very similar at the present day, but the global properties 
of the in situ halo differ substantially between the three galaxies in our study. In our 
simulations, the hierarchical buildup of structure is the driving force behind not only 
the accreted stellar halo, but also those halo stars formed in situ. 

Key words: methods: numerical - galaxies: format ion - galaxies: haloes- galaxies: 
structure 


1 INTRODUCTION 

Following Searle & Zinn (1978), much observational and the¬ 
oretical work on the Milky Way’s stellar halo has focused on 
the tidal stripping and disruption of satellite galaxies. The 
idea that galactic stellar haloes are built mainly by accre¬ 
tion is well supported by theoretical predictions of the stan¬ 
dard dark energy/cold DM (ACDM) cosmogony (White & 
Frenk 1991; Bullock & Johnston 2005; Cooper et al. 2010) 
and direct evidence of tidal streams around nearby galaxies 
(Belokurov et al. 2006; McConnachie et al. 2009; Martfnez- 
Delgado et al. 2010). However, recent work has shown that 
some aspects of the Milky Way’s stellar halo may be difficult 
to explain by accretion alone, notably its central concentra¬ 
tion and uniformity across the sky (Carollo et al. 2007; Bell 
et al. 2008; Cooper et al. 2011; ?; Deason et al. 2011; Helmi 
et al. 2011; Xue et al. 2011). 

Meanwhile, hydro dynamical simulations have predicted 
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a distinct Hn situ’ halo component, defined (loosely) as hav¬ 
ing formed bound to the Milky Way itself rather than to 
any of its hierarchical progenitors (Abadi et al. 2003; Brook 
et al. 2004; Zolotov et al. 2009; Font et al. 2011; Tissera 
et al. 2013; Pillepich et al. 2015). Such haloes are a natu¬ 
ral outcome of the ACDM model, which predicts that the 
vast majority of stars in a galaxy like the Milky Way form 
from the cooling of gas trapped by the galaxy’s own DM 
(DM) potential (White & Rees 1978; White & Frenk 1991). 
The bulk of these in situ stars can be identified with the 
kinematically cold, rapidly rotating Galactic disc, but the 
proto-Milky Way may also have suffered strong perturba¬ 
tions from satellites and quasi-secular rearrangement (e.g. 
‘disc flips’; Bett & Frenk 2012), or even wholesale destruc¬ 
tion and regrowth before the majority of present-day disc 
stars were formed (e.g. Sales et al. 2012; Aumer & White 
2013; Aumer et al. 2013). If real galaxies pass through such 
messy stages of formation, it seems likely that a significant 
fraction of stars formed in the early Galaxy would now have 
highly eccentric orbits. 
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Recent simulations find that these in situ processes cre¬ 
ate haloes that are more concentrated, metal-rich and oblate 
than those formed by accreted stars (McCarthy et al. 2012). 
This supports the hypothesis of a transition between an in 
situ and an accreted halo as an explanation for the appar¬ 
ently ‘bimodal’ properties of halo stars observed in the Milky 
Way (Carollo et al. 2010; Beers et al. 2012; Tissera et al. 
2014). 

However, quantitative results concerning the origin of 
in situ halo stars and their importance relative to accreted 
stars are still very uncertain. Where they rely on simula¬ 
tions, such conclusions can be particularly sensitive to the 
numerical methods used. Starting from identical initial con¬ 
ditions, current state-of-the-art simulations predict substan¬ 
tially different properties for the bulk of the in situ stellar 
mass in Milky Way-like DM haloes (e.g. Scannapieco et al. 
2012; Aumer et al. 2013), not just the few per cent that 
might be identified with an in situ halo. Moreover, the prop¬ 
erties of in situ haloes may be much more sensitive to certain 
modelling choices than those of massive stellar discs, includ¬ 
ing prescriptions for star formation and the treatment of the 
multi-phase interstellar medium (ISM; for example, the mix¬ 
ing of hot and cold in galactic winds, tidal streams and cold 
clumps of free-falling gas). 

Here we analyse the origin of in situ halo stars using 
three Milky Way-scale simulations run with the code de¬ 
scribed in Parry et al. (2012), one of the participants in the 
Aquila comparison project (Scannapieco et al. 2012). Two 
of the three DM haloes we simulate have also been simu¬ 
lated by Tissera et al. (2012, hearafter T12) and Tissera 
et al. (2013, 2014, hearafter T13, T14) using different ‘sub¬ 
grid’ recipes for star formation and feedback but an other¬ 
wise similar hydrodynamic solver and identical initial con¬ 
ditions. We define what we mean by in situ halo stars in 
a straightforward and easily reproducible way. Careful defi¬ 
nitions are particularly important for this problem because 
the concept of an in situ halo straddles an extremely fuzzy 
boundary between all the conventional Galactic components 
- disc, thick disc, bulge and halo. Based on these definitions, 
we discuss physical mechanisms by which in situ haloes are 
generated in our simulations, motivated by the fact that the 
mechanisms we identify are somewhat different from those 
previously discussed in the literature. In particular, we find 
that the growth of the in situ halo is very closely related to 
the accretion of satellites responsible for the growth of the 
canonical ‘ ex situ ’ halo. 

We proceed as follows. We describe our simulations in 
Section 2. In Section 3 we explain how we identify in situ 
halo stars and in Section 4 we examine their origins. Sec¬ 
tion 5 describes the present-day properties of our in situ 
halo. Section 6 investigates the satellite progenitor of in situ 
stars formed from stripped gas. In Section 7 we interpret our 
results, discuss the limitations of in situ halo models based 
on hydrodynamical simulations and compare with similar 
studies. A summary of our conclusions is given in Section 8. 
A detailed comparison with the results of T12 and T13 is in¬ 
cluded in Appendices A and B, along with a short discussion 
of numerical convergence. 


Table 1 . Numerical parameters adopted for the three simula¬ 
tions: DM and gas particle masses and the maximum gravita¬ 
tional softening length in physical units (defined as the scale of a 
Plummer kernel equivalent to the actual spline kernel used in the 
simulation). 



Mdm(Mq) 

Mgas (Mq ) 

e phys [P c ] 

Aq-C 

2.6 x 10 5 

5.8 x 10 4 

257 

Aq-D 

2.2 x 10 5 

4.8 x 10 4 

257 

Aq-E 

2.1 x 10 5 

4.7 x 10 4 

257 


2 SIMULATIONS 

We examine the stellar haloes that form in three smoothed 
particle hydrodynamics (SPH) simulations of Milky Way- 
mass galaxies. Dark-matter-only versions of these simula¬ 
tions formed part of the Aquarius project (Springel et al. 
2008) and we retain the Aquarius nomenclature for our three 
sets of initial conditions, labelling them Aq-C, Aq-D and 
Aq-E. The DM resolution (particle mass) in our SPH sim¬ 
ulations is similar to that of the ‘level 4’ simulation set in 
Aquarius. 

The Aquarius DM haloes were themselves extracted 
from a cosmological simulation in a cube of comoving vol¬ 
ume 100 3 Mpc 3 . They were chosen to have masses close to 
that of the Milky Way (~ 10 12 Mq) and to avoid dense en¬ 
vironments (no neighbour exceeding half the mass of the 
target halo within l/i -1 Mpc; Navarro et al. 2010). Ini¬ 
tial conditions for a resimulation of each halo were cre¬ 
ated with a ‘zoom’ technique, with higher mass bound¬ 
ary particles used to model the large scale potential and 
lower mass particles in an ~ 5/i -1 Mpc region surround¬ 
ing the target halo. Extra power was added to the initial 
particle distribution on small scales in the high resolution 
region, as described by Frenk et al. (1996). The numeri¬ 
cal parameters for each simulation, including the particle 
masses and gravitational softening lengths, are listed in Ta¬ 
ble 1. We assume a ACDM cosmology, with parameters 
Q m = 0.25, Qa = 0.75, Q h = 0.045, cr 8 = 0.9, n s = 1 
and Ho — 100/ikms -1 Mpc -1 — 73kms -1 Mpc -1 . 

Our simulation code is based on an early version of the 
PM-Tree-SPH code GADGET-3. Baryon processes are mod¬ 
elled as described in Okamoto et al. (2010) and Parry et al. 
(2012). Briefly, each gas particle represents an ISM with 
separate ‘hot’ and ‘cold’ (star-forming) phases. Gas parti¬ 
cles above a critical density (tih > 0.1cm -3 ) are assigned a 
cold phase mass according to their thermal energy and a lo¬ 
cal pressure according to a polytropic equation of state. Gas 
in the cold phase is converted to stars at a rate inversely 
proportional to a local dynamical time, which in turn de¬ 
pends on analytic approximations for the effective pressure 
and distribution of cold cloud sizes. Our code follows the nu- 
cleosynthetic production of individual elements separately, 
in particular the iron yields of Types II and I SNe. Mass, 
metals and energy returned by evolved stars and their SNe 
are smoothly distributed over 40 near neighbour gas parti¬ 
cles. The effective metallicity used to calculate the radiative 
cooling rate of a given particle is also a smoothed average. 
Kinetic energy is imparted directly to particles subject to 
SNe feedback, in proportion to the local velocity dispersion 
of DM. Particles are entrained in SN winds on a probabilis¬ 
tic basis and launched perpendicular to the plane of the 
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galactic disc. Wind particles are decoupled from the hydro- 
dynamical calculation until they reach an ambient density 
n h < 0.01 cm -3 or the time they have been decoupled ex¬ 
ceeds a limiting time. Parry et al. (2012) improve the treat¬ 
ment of fluid instabilities in the hydrodynamic scheme of 
Okamoto et al. (2010) by incorporating artificial conductiv¬ 
ity (Price 2008) and a time-step limiter (Saitoh & Makino 
2009); they also adjust the treatment of SN winds to re¬ 
duce the resolution dependence of mass loading and deposit 
thermal rather than kinetic energy from Type la SNe. The 
most relevant effect of these improvements is to suppress 
the formation of stars in cooling instabilities in diffuse cir- 
cumgalactic gas (e.g. Kaufmann et al. 2006; Keres et al. 
2012; Hobbs et al. 2013). In the results we present below, 
we explicitly identify stars formed through such instabilities 
and find that they contribute < 20 per cent of the mass 
of the in situ stellar halo in our simulations. This fraction 
depends somewhat on resolution (becoming larger at lower 
resolution) as we discuss in Appendix B. Other modelling 
uncertainties are also discussed in Section 7 and Appendix 
B. 

Details of two of our three hydrodynamical simulations 
(Aq-C and Aq-D) have been presented previously in studies 
that focused on satellite galaxies (Parry et al. 2012) and 
pseudo-bulge formation (Okamoto 2013). The luminosity 
function and luminosity-metallicity relation of the simulated 
satellites are comparable to those of dwarf galaxies in the 
Local Group. All three simulations result in galaxies with 
massive centrifugally supported discs as well as dispersion 
supported spheroids. 


3 SAMPLE DEFINITION 


The first step in defining our stellar halo sample is to iden¬ 
tify all stars belonging to the central (Milky Way-analogue) 
galaxy at the present day (redshift z = 0). We choose 
stars that lie within a radius 7*200 which encloses a sphere 
of mean density 200 times the critical value for closure 
( t 200 = 227 kpc for Aq-C and Aq-D, 202 kpc for Aq-E). 
From this sample, we isolate the halo by excluding stars 
that belong to satellite galaxies within 7*200 and stars that 
belong to the main galaxy disc or an inner spheroid (‘bulge’), 
as follows. 

Satellite DM haloes and their galaxies are isolated us¬ 
ing a version of the subfind algorithm (Springel et al. 2001) 
adapted by Dolag et al. (2009) to identify self-bound sub¬ 
structures, taking into account the internal energy of the gas 
when computing particle binding energies. All star particles 
bound to DM subhaloes at z = 0 are excluded from our halo 
star sample. 

The central galactic disc is identified by finding stars on 
orbits that are approximately circular and that lie close to 
a plane normal to the net angular momentum vector of the 
whole stellar component. A coordinate system is chosen such 
that the net angular momentum vector of all stars within 
0. 27*200 points in the positive ^-direction. The circularity of 
each star’s orbit is then defined as 


Jz 

Jcirc(E) ’ 


where J z is the z component of the star’s specific angular 
momentum and J c irc(E ) is the specific angular momentum 


Table 2. Total mass in the (outer) disc, ‘bulge’ and stellar halo 
regions of the (r, £ E ) plane according to our criteria. The final row 
gives the fraction of mass in the stellar halo region that is formed 
in situ. Our galaxies are roughly half the mass of the Milky Way; 
note that the disc mass quoted is only for stars with r > 5 kpc. 
Our ‘bulge’ definition includes all stars with r < 5 kpc, regardless 
of their kinematics (it does not correspond to a specific kinematic 
or photometric component in our simulations or in real galaxies). 


Aq-C Aq-D Aq-E 

Disc (r > 5 kpc, S E > 0-8) 3.0 4.3 4.3 

f109M i Bulge (r < 5 kpc) 34.5 21.9 25.2 

1 &) Halo (r > 5 kpc, S E < 0.8) 4.6 8.4 6.8 

In situ halo fraction 37% 33% 41% 


of a star with the same binding energy on a circular orbit. All 
stars with S E > 0.8 are identified with the central galactic 
disc and excluded from our halo star sample. 

Fig. 1 shows the distribution of stellar circularity as a 
function of radius in our three simulations. A concentra¬ 
tion of corotating stars on near-circular orbits extending to 
~ 30 kpc is obvious in all cases, which we identify with the 
thin disc. A number of streams on pro- and retrograde orbits 
are also visible at large radii. 

At r < 5 kpc the density of stars on non-circular orbits 
is comparable to the density of stars in the disc. We iden¬ 
tify this complex region with a galactic ‘bulge’ (using the 
term loosely, since this region is substantially more extended 
than the bulge of the Milky Way). To simplify our defini¬ 
tion of the stellar halo, we exclude all stars with r < 5 kpc, 
regardless of circularity. This cut is easy to apply to both 
models and data. It also follows the loose convention of most 
Milky Way stellar halo work, in which stars more than a few 
kiloparsecs interior to the solar neighbourhood are excluded 
(the exception being those high above the disc plane) even 
though the inward extrapolation of a canonical r~ 3 density 
profile would predict a substantial mass of halo stars in the 
centre of the Galaxy (see also the discussion in Cooper et al. 
2010 ). 

Fig. 1 further separates star particles into accreted 
(middle row) and in situ (bottom row) according to whether 
or not they are bound to the main branch progenitor of each 
DM halo at the first snapshot after their formation. Star par¬ 
ticles that are first bound to a DM halo other than the main 
progenitor are considered as accreted, even if they form in 
a subhalo of the main branch (i.e. if they form in a satellite 
galaxy of the Milky Way analogue) and are subsequently 
stripped 1 . Table 3 summarizes the total mass of the stellar 
halo and the relative proportion of in situ stars. 

Fig. 2 shows the mean density of accreted and in situ 
halo stars in spherical shells centred on the galaxy. The pro¬ 
file of the in situ halo has a similar shape and amplitude 
in all three simulations, with a slight steepening evident in 
the ‘bulge’ region of Aq-C. In both Aq-C and Aq-D, the 
accreted halo stars are less centrally concentrated than the 

1 This is an important difference with the work of T13, who in¬ 
cluded stars formed in bound satellites within 7*200 in the in situ 
halo as part of their ‘endo debris’ category. 
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r (kpc) r (kpc) r (kpc) 


Figure 1 . The distribution of stars in radius-circularity space for Aq-C (left column), Aq-D (centre column) and Aq-E (right column). 
Panels in the top row include all stars bound to the main DM halo (r < 90 kpc), while the middle and bottom rows include only accreted 
and in situ stars respectively. Dashed horizontal lines indicate the circularity cut used to define disc stars. Dashed vertical lines mark the 
5 kpc cut in radius used to define ‘bulge’ stars. All stars outside these regions are classified as halo stars. The colour scale corresponds 
to the logarithm of the number of star particles. 


in situ component, with a mild break due to accreted stars 
alone at 70 < r < 90 kpc, while in Aq-E the two components 
are almost indistinguishable. The accreted-m situ transition 
in these profiles at ~ 20 kpc is consistent with the average 
for Milky Way analogues in the GIMIC simulation (Font 
et al. 2011). 

Aq-D and Aq-E have ‘thick’ discs with a high degree of 
non-circular motion, apparent in the top row of Fig. 1 as a 
high density of stars at 0.5 < £e < 0.8 and 5 < r < 20 kpc. 
According to the aformentioned cuts on circularity and ra¬ 
dius, we classify these as halo stars. However, examining the 
variation of the circularity distribution with height above the 
disc plane reveals that these stars simply make up the low- 
circularity tail of a continuous distribution. The fraction of 
stars on circular orbits is highest close to the plane. 

The most important question from the point of view of 
this paper is not the origin of thick disc stars, but whether 
or not they can, or should, be distinguished as a separate 
galactic component. Our simulations differ greatly in this 


respect, in line with earlier studies of the origin and structure 
of thick discs in simulations (e.g. Sales et al. 2009). Fig. 1 
shows that accreted stars can make a significant contribution 
to the ‘disc’. In Aq-D, they contribute mainly to the ‘thin’ 
disc - the thick disc is formed in situ. In Aq-E, accreted and 
in situ ‘disc’ stars contribute at a similar ratio over a wide 
range of circularity and radius. 

In the context of observations of the Milky Way, a ge¬ 
ometrically oblate stellar component, intermediate between 
disc and halo, was identified by Yoshii (1982) and Gilmore & 
Reid (1983), but optimal and objective ways to classify this 
component are still under debate (see, for example, recent 
reviews by Ivezic et al. 2012; Rix &; Bovy 2013). Classifi¬ 
cations have been suggested based on metallicity and kine¬ 
matics that imply different distinctions between thick disc 
and halo stars (Ivezic et al. 2008; Bovy et al. 2012; Rix & 
Bovy 2013; Schlesinger et al. 2014; Ruchti et al. 2015). The 
chemodynamical distribution functions of our three galaxies 
are different from one another, and almost certainly differ- 
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logio r /kpc 

Figure 2. Spherically averaged stellar density profiles for Aq- 
C (top), Aq-D (centre) and Aq-E (bottom). The red dashed and 
green dot-dashed lines correspond to the accreted and in situ halo 
components respectively. The solid black line is the sum of these 
and can be compared to measurements for the Milky Way (black 
star; Fuchs & Jahreifi 1998; Gould et al. 1998). The blue dotted 
line corresponds to disc stars, selected by their orbital circularity. 
The minimum extent of the radial axis is set by the gravitational 
softening length for star particles. A grey dashed vertical line at 
5 kpc marks the ‘bulge’ region that we exclude from our analysis 
of the stellar halo. 


ent from the Milky Way distribution function on which these 
observational definitions are based. Therefore, we do not be¬ 
lieve that additional cuts to separate a thick disc component 
from the thin disc and halo (for example, using a threshold 
metallicity, a limiting rotation velocity or vertical extent) are 
helpful for the interpretation of our simulations. This should 


Table 3. Breakdown of all in situ halo stars into three subtypes, 
according to their formation mechanism. 



Aq-C 

Aq-D 

Aq-E 

Heated disc 

2 .8% 

26.0% 

31.0% 

Stripped gas 

59.8% 

56.7% 

56.9% 

Smoothly accreted gas 

37.3% 

17.3% 

12 .1% 


be kept in mind when comparing our results to observations 
of the Milky Way. 


4 THE ORIGIN OF IN SITU STARS 

In this section we look in more depth at the origin of the in 
situ component of the stellar halo. The density of this com¬ 
ponent exceeds that of accreted halo stars in the inner ~ 20 
kpc of our galaxies. It may thus be very important for spec¬ 
troscopic observations of halo stars in the solar neighbour¬ 
hood and in surveys of main-sequence turnoff stars within a 
few kiloparsecs of the Milky Way disc plane. 

In order to trace how in situ stars formed in our simu¬ 
lations, we separate them into three disjoint subcategories: 

(i) ‘heated disc’ stars, which met the thin disc circularity 
criterion when they were formed, but are not in the disc at 
z — 0; 

(ii) stars formed from ‘stripped gas’, brought into the 
main DM halo bound to a subhalo and subsequently 
stripped by tidal forces or ram pressure; 

(iii) stars from ‘smoothly accreted gas’, which enters the 
main DM halo through direct (smooth) accretion. 

Categories (ii) and (iii) are easily distinguished by trac¬ 
ing the DM halo membership history of the parent gas par¬ 
ticle for each star particle. The gas from which heated disc 
stars form must originally have either been stripped from 
a subhalo or smoothly accreted, so these stars could also 
be classified into the second or third categories. However, 
in this case it is the fact that they formed in the thin disc 
and were scattered out of it, rather than how their parent 
gas particle arrived in the disc, that we consider to be most 
important. 

The fraction of in situ halo stars in each category is 
shown in Table 3. It is clear that there is a large variation 
between the three simulations, although the in situ stars 
forming from gas stripped from satellites dominate in all 
cases. In the next section we discuss each category in more 
detail and compare their properties. 

4.1 Heated disc stars 

The central galaxies in our three simulations undergo several 
episodes of disc destruction and regrowth at z > 3. Over the 
redshift range 3 > z > 2, a stable disc is established. This 
disc continues to grow until z = 0, although its angular 
momentum axis may precess. Our heated disc category only 
includes stars that once belonged to this stable disc. Stars 
on highly circular orbits may be scattered to more eccentric 
orbits by secular evolution and satellite impacts (e.g. Purcell 
et al. 2010). We refer to this loosely as ‘heating’, in the 
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Table 4. Properties of the present-day thin stellar discs in our 
three simulations, defined by circularity £e > 0.8. Columns from 
left to right give formation redshift, mass fraction in place at 
formation redshift, and redshift at which half the z = 0 mass is 
in place. 



*form 

/form 

*1/2 

Aq-C 

2.54 

12% 

0.92 

Aq-D 

2.32 

8% 

0.83 

Aq-E 

2.20 

16% 

0.76 


sense of an increase in non-circular motion. These perturbed 
disc stars are likely to have a clear kinematic and chemical 
relationship to those in the present-day thin disc. 

We identify all star particles in the z — 0 disc, as defined 
in Section 3, that exist in a given earlier snapshot and use 
these to define J z (assuming that the number of star parti¬ 
cles scattered into the disc is negligible). We then apply the 
circularity threshold £e > 0.8 to identify all newly formed 
star particles in the disc at that snapshot. Any of these that 
have £e < 0.8 at z = 0 are assigned to our heated disc 
category. 

Beyond a certain redshift, Zf or m, we can no longer reli¬ 
ably identify a stable progenitor of the z = 0 disc, and hence 
we cannot define J z . This limit is due to the small number 
(< 100) of ancient disc stars and increasing frequency of fluc¬ 
tuations in the central potential that destabilize protodiscs. 
Table 4 gives Zf orm for each of our simulations, along with 
Z\/ 2 , the redshift by which the z = 0 disc has assembled half 
its final mass. Table 4 also gives /form, the mass fraction of 
z — 0 disc stars that form earlier than Zf or m • This fraction is 
no more than 16 per cent (Aq-E). Stars scattered from this 
unidentified protodisc (and any others discs that were com¬ 
pletely destroyed before Zf or m) are considered to fall into one 
of the other two in situ categories, according to the origin 
of their parent gas particle. 


4.2 Stars from stripped gas and smoothly 
accreted gas 

Halo stars can form directly in the circumgalactic medium, 
either in quasi-free-falling cold gas clouds (not associated 
with DM clumps) or the gaseous tidal or ram pressure 
stripped streams of satellite galaxies. We distinguish be¬ 
tween these two possibilities based on whether or not the 
parent gas particle of a given star particle was bound to an¬ 
other DM halo before being bound to the main halo. Stars 
forming from stripped satellite gas particles may be chemi¬ 
cally and kinematically similar to stars in the accreted stellar 
halo. In contrast, stars forming in gas condensing out of the 
hot hydrostatic gas halo, or other ‘smoothly’ accreted cold 
clumps, may have properties more similar to those expected 
of an in situ halo formed by monolithic collapse. 

Fig. 3 shows the absolute star formation rate of each 
in situ category as a function of time elapsed since the big 
bang. These star formation rates are low compared to those 
typical of the stable disc and the progenitors of accreted 
stars (~ lMQyr -1 ). The majority of halo stars that form 
in stripped or smoothly accreted gas are more than 9 Gyr 
old, only marginally younger than the typical age of accreted 



Figure 3. Formation history for in situ halo stars assigned to 
each of our three in situ formation mechanisms. Time is measured 
from the origin of the universe. Peaks in the formation of stars in 
smoothly accreted gas are clearly correlated with those in stripped 
gas at t < 4 Gyr. Vertical dotted lines mark 2 :f orm and 2q/ 2 as 
given in Table 4. 

stars. In Aq-D and Aq-E, there are also ~ 2 Gyr-long bursts 
of in situ star formation starting at t « 4 and 6 Gyr respec¬ 
tively. Interestingly, these also correspond to episodes of for¬ 
mation for scattered disc stars (blue) and thin disc stars (not 
shown). These episodes correspond to the rapid infall of cold 
gas on to the disc during periods in which several relatively 
massive satellites are being disrupted simultaneously. 

Another notable feature of Fig. 3 is that the star forma¬ 
tion rate in smoothly accreted gas (red) is clearly correlated 
with that in stripped gas (green), especially in the domi¬ 
nant early epoch of in situ halo formation (ages > 9 Gyr). 
This correlation persists even if we select only stars forming 
at r > 30 kpc, far away from the disc, suggesting that the 
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conditions under which most in situ halo stars form are in 
fact related to the accretion and stripping of gas-rich satel¬ 
lites. It appears that star formation may be triggered by the 
mixing of free-floating gas from the hydrostatic halo with 
star-forming stripped gas. We also see corresponding peaks 
in the accreted halo star formation rate, suggesting that star 
formation is triggered in the infalling satellites as well. 

These effects must arise from the subgrid model of star 
formation. They may therefore occur, to a greater or lesser 
extent, with other similar star formation prescriptions used 
in the literature. Since there are no unambiguous observa¬ 
tional tests of the predictions of models for star formation 
in the diffuse circumgalactic gas, it is very difficult to judge 
whether or not the behaviour of a particular simulation is 
physically plausible (see also the discussion in Appendix B). 
Isolating dependencies on modelling uncertainties such as 
these is one of our motivations for separating simulated halo 
stars into categories based on their physical origin. 

Regardless of how they actually form, the classification 
of gas particles as ‘smoothly’ accreted is resolution depen¬ 
dent: gas bound to a low-mass infalling DM halo, and so 
classified as ‘stripped’ at high resolution, would be classified 
as ‘smoothly accreted’ in a lower resolution simulation that 
does not resolve the parent halo. In Appendix B we con¬ 
clude that, in our simulations, the uncertainty in the total 
mass of the smoothly accreted halo component is dominated 
by this classification uncertainty, rather than the resolution 
dependence of star formation efficiency. 


5 IN SITU HALOES AT z = 0 

In this section, we examine the observable characteristics of 
in situ halo stars at the present day, starting with a summary 
of halo properties and then looking in more detail at regions 
analogous to the solar neighbourhood. 

5.1 Whole halo 

Fig. 4 compares the spherically averaged density profiles of 
our three in situ halo categories and accreted halo stars. 
In the r < 20 kpc region where in situ halo stars dominate 
over accreted stars, they contribute roughly equal mass frac¬ 
tions; the exact proportions vary from halo to halo. We see a 
strong correspondence between stars formed from stripped 
and smoothly accreted gas at all radii, which, in combination 
with Fig. 3, suggests that they form with a similar distribu¬ 
tion in both space and time. As expected, heated disc stars 
have a steeper profile, with most concentrated at r < 20 kpc. 

Fig. 5 shows Toomre diagrams (Sandage & Fouts 1987) 
that compare the amplitude of circular and radial motion 
for different components. A galactocentric UVW velocity 
frame (e.g. Binney & Merrifield 1998, p. 627) is defined with 
respect to the thin disc in each simulation. 

Of the three simulations, the stellar halo in Aq-C has 
kinematic properties most similar to those measured for 
the Milky Way. The peak rotational velocity of the disc 
is ~ 220kms _1 . The heated disc stars (blue) rotate in 
the same sense as those on circular orbits, with a lag of 
~ 40-180 km s -1 . Stars formed from stripped and smoothly 
accreted gas are kinematically indistinguishable from each 
other, once again pointing to a close correlation between the 
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Figure 4. Spherically averaged density profiles of halo stars in 
the accreted and three in situ components for Aq-C (top), Aq-D 
(centre) and Aq-E (bottom). 

dynamics of the two components. In the Toomre diagram 
they resemble the classical Milky Way halo, with zero net 
rotation and high radial velocity dispersion. Accreted halo 
stars show a similar distribution overall, but with notable 
overdensities due to individual streams, some of which have 
a net retrograde motion. 

The heated disc stars in Aq-D and Aq-E have similar 
kinematics to those in Aq-C, but the stripped/smooth in 
situ haloes have a greater net rotation. In Aq-E, all three 
components once again resemble one another, although the 
stripped- and smooth-gas halo stars have a greater velocity 
dispersion. An underlying stripped/smooth in situ halo may 
still be present, but the bulk of in situ halo stars are more 
similar kinematically to the Milky Way thick disc. The be- 
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Figure 5. Toomre diagrams of the whole stellar halo. For 
stripped-gas, smooth-gas and heated disc in situ halo stars (red, 
green and dark blue respectively), contours mark the regions en¬ 
closing 10, 30, 50, 70 and 90 per cent of the stellar mass. For the 
accreted halo (cyan) only 10, 50 and 90 per cent levels are shown. 
The dashed vertical line marks the rotation velocity of the disc 
at 8 kpc. 



haviour of accreted stars once again resembles that of the in 
situ component, even to the extent that they have a strong 
prograde rotation in Aq-E. Accreted halo components with 
prograde rotation were noted by Abadi et al. (2003) and also 
found in Milky Way-like systems in the GIMIC simulations 
(Font et al. 2011; McCarthy et al. 2012) 

Finally, in Fig. 6, we examine the normalized metallicity 
distribution functions (MDFs) of each component of the in 
situ halo. Heated disc stars have the highest median [Fe/H] 
and narrowest dispersion. Their MDF resembles that of the 
thin disc, but is slightly more metal poor (by ~ 0.5 dex in 


Aq-C). Both in situ and accreted halo stars are systemati¬ 
cally more metal poor than heated disc stars. 

The MDF of in situ halo stars formed from stripped gas 
is very similar to that of accreted satellite stars, with a me¬ 
dian systematically higher by no more than 0.1 dex. This is 
to be expected, as the dense cold gas stripped from satellites 
will have been enriched by the same stellar populations that 
make up the accreted halo. Moreover, very similar distribu¬ 
tions will also result if prolonged star formation occurs in 
satellite galaxies while their gas is being stripped. The over¬ 
all in situ MDF is close to that of the stripped-gas stars, 
since they dominate the in situ mass budget. 

Looking in detail, the degree of similarity between the 
MDFs of the various components varies in each of our three 
simulations. This may depend on the extent to which the 
satellite galaxies contributing the bulk of stripped-gas stars 
are the same as those that contribute the majority of ac¬ 
creted stars. Since gas can be more easily expelled from 
shallower potentials, the most massive and metal-rich ac¬ 
creted progenitor galaxies are likely to retain the most gas 
when they enter the main DM halo. Stars stripped from 
these galaxies are expected to dominate the accreted halo, 
particularly near the centre. We investigate the relative con¬ 
tributions of star-forming stripped gas and directly accreted 
stars from different progenitors in Section 6 below. 

Of the different in situ components, it is the stars that 
formed from smoothly accreted gas that have the lowest me¬ 
dian metallicity and the broadest dispersion. This is con¬ 
sistent with the expectation that the gas surrounding each 
galaxy will be a mix of its own metal-rich ejecta and a large 
quantity of ‘pristine’, or only marginally enriched, gas from 
direct cosmological infall (Crain et al. 2010). In our simula¬ 
tions, the MDF is the only clear distinction between ‘smooth 
gas’ stars and ‘stripped gas’ halo stars. 


5.2 The Solar Neighbourhood 

As a rough analogue of the solar neighbourhood region most 
relevant to current observations, we examine the average 
properties of halo stars in a torus of cross-sectional diameter 
4 kpc and galactocentric radius r — 8 kpc in the plane of 
the thin disc. 

Table 5 summarizes the fraction of stars in each compo¬ 
nent. For a more direct comparison to the real data, we have 
grouped heated disc stars and stars that meet our thin disc 
circularity cut into a single disc component, because the typ¬ 
ically high circular velocities of heated disc stars would most 
likely result in them being classified as ‘thick disc’ rather 
than halo stars in observations. Approximately 10 per cent 
of the stellar mass then remains in a component resembling 
the ‘classic’ halo, of which accreted stars contribute between 
34 and 67 per cent. 

Toomre diagrams in this region are almost identical to 
those in Fig. 5. The biggest differences in comparison to the 
overall halo are found in the solar neighbourhood MDFs, 
which are shown in Fig. 7. Stars formed from smoothly ac¬ 
creted gas that end up in the solar neighbourhood are more 
metal rich on average, such that their MDF has a very sim¬ 
ilar shape and amplitude to the accreted halo. This may be 
because the metal-poor contribution of this component seen 
in Fig. 6 is dominated by stars forming at large radii from 
gas that has not been polluted by the galactic wind of the 
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Table 5. Breakdown of all stars in the solar neighbourhood. Heated disc stars are grouped together with thin disc stars in this table. 
The top two rows give fractions of total stellar mass, while the lower three rows give fractions of stellar halo mass (second row) only. 


Mass (10 s M 0 ) 

Aq-C 

Aq-D 

Aq-E 

Disc stars (thin + thick) 

26.4 (92.3%) 

27.4 (88.7%) 

28.0 (83.0%) 

Halo stars 

2.21 (7.7%) 

3.49 (11.3%) 

5.72 (17.0%) 

Accreted 

0.747 (33.9%) 

1.79 (51.1%) 

3.83 (67.0%) 

Stripped gas 

0.837 (37.9%) 

1.23 (35.3%) 

1.53 (26.7%) 

Smoothly accreted gas 

0.618 (28.0%) 

0.476 (13.6%) 

0.356 (6.23%) 



[Fe/H] 

Figure 6. [Fe/H] distributions for the disc, accreted halo and 
three in situ halo components. Distributions are normalized by 
the total mass of stars in each component. 



Figure 7. MDFs, as Fig. 6, but here in the solar neighbourhood. 
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central galaxy. Other components have the same relation¬ 
ship to one another as those in Fig. 6. Hence, we find no 
substantial differences between the properties of the in situ 
halo in the solar neighbourhood and the in situ halo overall. 
This is not surprising because we have already seen that the 
bulk of the in situ halo is concentrated within r < 20 kpc. 


6 SATELLITE PROGENITORS 

We have shown that our simulated stellar haloes are dom¬ 
inated by satellite accretion: the bulk of halo stars are 
stripped directly from satellites, and the majority of l in situ? 
halo stars form from gas stripped from satellites. However, 
as Fig. 2 demonstrates, stars formed in situ from stripped 
gas have a more centrally concentrated spatial distribution 
at z = 0 than directly accreted stars. In this section we ex¬ 
amine the satellites which contribute to the halo, with the 
aim of determining how their infall times, masses and bary- 
onic content affect the spatial distribution of the in situ and 
accreted components. 

We first ask whether the subset of satellite progenitors 
contributing the gas from which an in situ stellar halo forms 
is the same subset contributing accreted stars. Fig. 8 com¬ 
pares the mass fractions of stars formed from stripped gas 
and directly accreted stars associated with each progenitor 
satellite. In all three simulations, satellites that contribute 
significantly to one component also tend to contribute signif¬ 
icantly to the other. There is substantial variation in detail 
between the three haloes, reflecting their different accretion 
histories. A larger scatter is apparent in Aq-D, as well as a 
noticeable fraction of gas-poor contributors (lower-right area 
of the plot) relative to Aq-C and Aq-E. A larger fraction of 
those satellites also survive to z = 0 without being disrupted 
(blue points). The smaller number of surviving satellites in 
Aq-C reflects a quieter recent merger history. 

In Fig. 9 we isolate the top three satellite progenitors of 
the stripped-gas in situ component and plot the density pro¬ 
files of the accreted and stripped-gas stars they contribute. 
Aq-E stands out as having the most similar profiles for the 
two components, both of which are slightly steeper than the 
total accreted profile. In this case, the satellites plotted ac¬ 
count for around 40 per cent of the total stripped-gas halo. 
With the possible exception of accreted stars in Aq-C, both 
accreted and stripped-gas stars from the top three satellites 
are distributed like the bulk of the stellar halo. 

Fig. 9 suggests that the greater central concentration of 
the stripped-gas halo profile is not simply because most of 
the progenitor gas particles originate in more massive pro¬ 
genitors, which sink more rapidly through the action of dy¬ 
namical friction. If that were the case, we might also expect 
stars accreted from the same progenitors to be more cen¬ 
trally concentrated than the accreted halo overall. Instead, 
we see that the debris profiles of the three most massive in¬ 
dividual progenitors are very similar to the profile of the en¬ 
tire accreted halo. In our simulations, at the present day, ac¬ 
creted stars are distributed (on average) over the same range 
of radii at which they were liberated from their parent satel¬ 
lites. Conversely, we have confirmed that the stripped-gas 
particles from which in situ halo stars form dissipate some 
of their orbital energy between the times of stripping and 
star formation. The present-day distribution of the stripped 
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Figure 8. Mass fraction of the accreted halo and stripped-gas 
halo contributed by disrupted (red) and surviving (blue) satellite 
galaxies. The sizes of the points are proportional to the logarithm 
of the satellite’s total mass at infall, as shown by the legend in the 
first panel. The diagonal dashed line indicates an equal fractional 
contribution to the accreted and stripped-gas components. 


gas halo is therefore imprinted at the time those stars form, 
rather than the time at which their parent gas particles are 
stripped. 


7 DISCUSSION 

Our finding of accretion-triggered star formation in 
smoothly accreted gas may be, at least in part, a conse¬ 
quence of our hydrodynamics scheme. In some hydrody¬ 
namic models, cold clouds can condense independently of 
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Figure 9. Spherically averaged density profiles of halo stars asso¬ 
ciated with the top three contributors to the ‘stripped-gas’ com¬ 
ponent. The profiles of stars formed from gas contributed by those 
galaxies are shown as dashed, cyan lines, while the profile of stars 
accreted from those satellites directly are shown as solid, green 
lines. Thinner, black lines of the same styles show the correspond¬ 
ing total profiles. 


tidally stripped gas, while in others, condensation and mix¬ 
ing by fluid instabilities may both be suppressed (e.g. Kauf- 
mann et al. 2006; Keres et al. 2012; Parry et al. 2012; Hobbs 
et al. 2013). As we noted in Section 4.2, there are no strong 
observational test of different models at this level of detail. 
In Appendix B, we carry out a simple resolution study of 
the smoothly accreted gas component. At lower resolution, 
more baryons are retained by the main halo and relatively 
more of these baryons are classified as smoothly accreted. 


The overall efficiency with which smoothly accreted gas is 
converted to stars also varies systematically with resolution, 
although (in our simulations) this is a minor contribution to 
the change of in situ stellar halo mass, relative to the first 
two effects. 

Given such uncertainties, it is worth noting that, with 
a single SPH code, the different cosmological initial condi¬ 
tions of our three simulations result in markedly different in 
situ halo properties at z = 0. It is obvious that the mass 
spectrum, arrival time distribution and prior star formation 
histories of accreted satellites directly determine the proper¬ 
ties of the accreted halo. Remarkably, however, in our SPH 
simulations, these factors are also extremely important in 
the formation of the in situ halo. Without satellite accre¬ 
tion, there would be hardly any in situ halo stars in our 
simulations. Moreover, it is unlikely that all the stars we 
classify as heated disc halo stars would be considered as 
such by observers. The vast majority are on orbits of mod¬ 
erate circular velocity, corotating and coplanar with the thin 
disc, and would likely be classified as ‘thick disc’ stars as dis¬ 
cussed in Section 3. In Aq-C, a large fraction of these stars 
belong to the end of a bar in the galactic plane that extends 
slightly beyond our ‘bulge’ cut at 5 kpc. 

Our model is clearly not unique. Nonetheless, our in situ 
haloes are compatible with the basic observable properties 
of the Milky Way’s stellar halo, if we exclude scattered disc 
stars. Stars formed from stripped and smoothly accreted gas 
in the main halo are old and metal poor, with a significant 
tail to low metallicity. The kinematic and chemical similar¬ 
ities are closest in Aq-C, which also has the most Milky 
Way-like thin disc. 

Two of our haloes conform to the expectation from ear¬ 
lier work that in situ halo stars are more concentrated than 
accreted stars, but in the third, the in situ and accreted halo 
density profiles are almost identical. Interestingly, it is only 
in this third case that there is a significant difference in the 
MDFs of in situ and accreted stars. 

McCarthy et al. (2012) studied the origin of in situ halo 
stars in 412 Milky Way-mass galaxies from the GIMIC sim¬ 
ulations, which have a DM particle mass ~ 300 times larger 
than our simulations. With the same definition of disc stars 
and accreted stars, they find in situ halo mass fractions rang¬ 
ing from 20 to 60 per cent. Their in situ stars are, on aver¬ 
age, younger and more metal rich than their accreted halo 
stars, with a centrally concentrated oblate distribution and 
prograde rotation. They demonstrate that these properties 
result from star formation at z ~ 1 in ‘proto-discs’ that are 
subsequently destroyed. 

Although our three simulations produce total in situ 
mass fractions consistent with the distribution found by Mc¬ 
Carthy et al. (2012), there are some notable differences in 
the origin of in situ stars. Accreted gas is the dominant con¬ 
tributor to our in situ haloes, whereas McCarthy et al. find 
that approximately half of their in situ halo mass is formed 
from gas that is shock-heated to the virial temperature of 
the main DM halo. Our in situ DM haloes form at some¬ 
what higher redshift and we can trace back our z = 0 discs 
to z ~ 2.5, supporting the suggestion of McCarthy et al. 
that the Aquarius simulations have quieter-than-average ac¬ 
cretion histories. As McCarthy et al. do not compute the 
fraction of gas forming their in situ component that was 
previously bound to satellites, it is likely that their proto- 
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disc stars mix together both stripped and smoothly accreted 
gas, and possibly some fraction of our heated disc popula¬ 
tion. Font et al. (2011) note that formation from stripped 
gas is apparent in the same set of galaxies from the GIMIC 
simulations and is likely responsible for their in situ stars 
at r > 20 kpc. We find that this mode of formation domi¬ 
nates at somewhat smaller radii. This may be a consequence 
of the better spatial resolution in our simulations; a smaller 
gravitational softening length means satellite gas tends to be 
more tightly bound, allowing it to sink further into the po¬ 
tential before being stripped by ram pressure or tidal forces 
(see also Parry et al. 2012). 

As noted in the introduction, a series of papers by Tis- 
sera et al. (T12, T13, T14) examine in situ halo stars in 
six SPH simulations from the Aquila suite (Scannapieco 
et al. 2012). They also use gadget-3 for these simulations, 
but with different implementations of subgrid physics and 
slightly lower resolution (a factor of 5 in DM particle mass; 
Scannapieco et al. 2012). Two of their simulations, Aq-C 
and Aq-D, have the same initial conditions as the simula¬ 
tions with the corresponding labels in this paper. To aid the 
reader in comparing our results with this series of papers, 
we comment here on the differences between our sample def¬ 
initions and those of Tissera et al. 

T12 separated disc and halo stars with a circularity cut, 
Jz/Jcirc(E) > 0.65, which classifies many more stars as ‘disc’ 
than our cut. Also, they did not restrict the height of disc 
stars above the mid-plane. They excluded bulge stars with 
an energy criterion that roughly equates to a radius of 5 kpc, 
comparable to our ‘bulge’ cut (see T12). 

In the analysis of Tissera et al., in situ halo stars are 
divided into ‘inner’ (more bound) and ‘outer’ (less bound) 
populations according to a cut in relative binding energy. 
This equates to a radial cut in the range 14 < r < 36 kpc 
(15 and 19 kpc for Aq-C and Aq-D respectively). They fur¬ 
ther subdivided these populations by origin, according to 
whether the stars were formed outside the virial radius (‘de¬ 
bris’) or inside (‘endo-debris’). This definition is different 
from ours, which counts a present-day halo star particle as 
‘accreted’ if it formed bound to a satellite subhalo, regard¬ 
less of whether the satellite was inside or outside the virial 
radius of the main DM halo at the time. This may not be sig¬ 
nificant, because the fraction of accreted halo stars formed 
in bound satellites after their infall is typically small (see 
McCarthy et al. 2012). 

The definition of heated disc halo stars in T13 requires 
the star particle to have J z /Jdrc(E) > 0.65 at its formation. 
Our definition also requires those stars to have been scat¬ 
tered from the disc that survives at z = 0; we would classify 
some fraction of their heated disc stars as formed from ei¬ 
ther stripped or smoothly accreted gas. This may explain 
why T13 find a much larger heated disc fraction in Aq-C 
despite their less stringent disc circularity cut (indeed, they 
note that heated disc stars have retrograde rotation in one 
of their simulations). 

Allowing for these differences of definition between our 
study and that of T13, a number of similarities are clear. The 
dominance of in situ over accreted halo stars in our simu¬ 
lations at r < 20 kpc is also seen in their ‘inner halo’ pop¬ 
ulations. Moreover, they also find that scattered disc stars 
make a negligible contribution outside this region. The kine¬ 
matic properties, as shown in our Fig. 5, are similar. In Aq-C 


and Aq-D, T13 report a median [Fe/H] of endo-debris stars 
lower than those of debris stars by 0.17 and 0.26 dex, re¬ 
spectively. This is the opposite of what we find, but could 
be easily explained if the T13 endo-debris definition includes 
a similar or larger fraction of stars formed from smoothly ac¬ 
creted gas. This is supported, at least qualitatively, by our 
findings in Appendix A, which indicate first that the star 
formation efficiency for the central galaxy is much higher in 
the Tissera et al. simulations, and secondly that stars in our 
simulations meeting the Tissera et al. endo-debris definition 
are predominantly formed in satellites after their accretion 
by the main halo. 

Recently, Pillepich et al. (2015) described the in situ 
stellar halo of the Eris simulation (a zoom of an isolated 
halo, M 200 = 8x 10 11 Mq, forming a late-type galaxy with 
M* = 3.9 x 10 10 Mq; Guedes et al. 2011). The particle mass 
of Eris is ~ 1/3 of that in our simulations. They find a 
transition in the density profile between in situ stars (in¬ 
cluding the disc) and accreted stars at ~ 10 kpc. Having 
cut out a cylinder of height ±5 kpc and radius 15 kpc to ex¬ 
cise the disc and bulge, they define an inner halo (all stars 
r < 20 kpc excluding this cylinder) and an outer halo (all 
stars 20 < r < 235kpc). They find that ~ 25 per cent of 
stars in their inner halo were formed in situ , falling to only 
3 per cent in their outer halo. Applying the same cuts to our 
haloes C, D and E, we find in situ fractions of (83, 78, 76) per 
cent for an equivalent inner halo region and (29,31,65) per 
cent for an equivalent outer halo region. These much larger 
fractions are most likely because Pillepich et al. count both 
disc and halo stars in their definition of in situ stellar mass 
- our discs are larger than the disc of Eris and hence are not 
excised completely by their geometric cut. If we count only 
our halo stars in these regions (excluding disc stars by their 
circularity as above), the in situ mass fractions for our three 
simulations are (20,18,15) percent in the Pillepich et al. in¬ 
ner halo and (14,11,17) per cent in the outer halo. Larger 
‘outer halo’ in situ mass fractions in our simulations may 
well result from differences between our subgrid models and 
those used for Eris , in particular their higher density thresh¬ 
old for star formation and neglect of metal line cooling in 
gas with temperature > 10 4 K (Guedes et al. 2011). 


8 CONCLUSIONS 

We find the following properties of in situ stellar halo stars 
in three high resolution SPH simulations of Milky Way ana¬ 
logues: 

(i) The in situ stellar halo accounts for 30-40 per cent of 
the stellar mass outside the thin disc (£e > 0 .8) and inner 
spheroid (‘bulge’; r > 5 kpc). This fraction includes stars 
that observers may classify as belonging to a thick disc. 

(ii) The in situ halo dominates over the accreted stellar 
halo at r < 30 kpc in two out of our three simulations. In 
the third simulation, both in situ and accreted halo stars 
have almost the same volume density distribution. 

(iii) Between 2 and 30 per cent of the in situ halo com¬ 
prises stars scattered from near-circular orbits in the plane 
of the thin disc to more eccentric orbits. These form at 
lookback times of 5 to 9 Gyr, retain high circular veloci¬ 
ties and have relatively narrow MDFs with median — 1 < 
[Fe/H] < —0.5. This component resembles the Milky Way’s 
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thick disc, although we note that, in two out of our three 
simulations, kinematically selected eccentric/thick discs also 
have a significant contribution from both accreted stars and 
stars formed from accreted gas. 

(iv) The rest of the in situ halo stars form in gas clouds in 
the circumgalactic halo, on highly eccentric orbits near the 
disc plane, or in the chaotic gas-rich stage of the galaxy’s 
formation (z > 2.5), before a stable disc is established. In 
one simulation, these halo stars strongly resemble the ‘classi¬ 
cal’ isotropic, metal-poor, dispersion-supported stellar halo 
of the Milky Way. In the other two simulations, most such 
halo stars have significant rotation in the same sense as the 
thin disc. 

(v) We identify two distinct origins for these stars: gas 
that has been stripped out of satellite galaxies by tides and 
ram pressure, and gas that is incorporated directly into the 
‘smooth’ halo of the main galaxy by cosmological infall and 
SN-driven outflow from the central galaxy. 

(vi) Stars formed in situ from stripped gas have a very 
similar MDF to the accreted stellar halo, because this gas 
is brought in by the same progenitor satellites. Halo stars 
formed from smoothly accreted gas have a broader MDF 
and are, on average, the most metal poor of the components 
we identify. 

(vii) The density profile of halo stars formed in situ from 
stripped gas is more concentrated than that of the stars ac¬ 
creted from the same progenitor galaxies. This reflects dissi¬ 
pative collapse of this stripped gas after it is liberated, rather 
than differences in how and when star and gas particles are 
stripped, or in the contributions of different satellites. 

(viii) In all cases, among the in situ halo stars not scat¬ 
tered from the disc, there is almost no difference in the 
present-day phase-space distribution of those formed in 
stripped and smoothly accreted gas. The correspondence is 
so close that we suggest that, in our simulations, star for¬ 
mation in the hot gaseous halo is directly triggered by the 
passage of dense clumps of star-forming stripped gas. 

(ix) The properties of in situ stars in the solar neighbour¬ 
hood are representative of the in situ halo overall, except 
that stars formed from smoothly accreted gas in this region 
are notably more metal rich. 

Based on these findings, we conclude that essentially all 
halo stars in our simulation are the result of cosmological 
accretion and merging, with no obvious bimodality due to 
an in situ halo forming through a quasi-monolithic collapse 
of enriched diffuse gas. 

Comparison with other studies that describe in situ 
stellar haloes in hydrodynamical models of Milky Way-like 
galaxies suggests that a wide variety of physical mechanisms 
can be responsible for their formation. The relative impor¬ 
tance of these processes is particularly sensitive to the im¬ 
plementation and resolution dependence of subgrid physi¬ 
cal prescriptions for star formation and feedback. An im¬ 
proved understanding of these effects and observational con¬ 
straints on star formation in free-floating gas clouds are thus 
necessary before definitive conclusions about in situ stellar 
halo formation can be drawn from any particular simulation. 
This uncertainty is currently the most important limitation 
on theoretical predictions for the in situ stellar halo. The 
comparison of simulations from different groups will be ex¬ 


tremely useful in this regard, hence our emphasis on simple 
origin-based definitions for the different halo components. 
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APPENDIX A: QUANTITATIVE COMPARISON 
WITH TISSERA ET AL. 

The initial conditions of haloes Aq-C and Aq-D have already 
been simulated by another group using different subgrid pre¬ 
scriptions and at slightly lower resolution, as reported by 
Tissera et al. (2012, 2013, 2014). This appendix provides a 
more detailed quantitative comparison of results for these 
two haloes to supplement the discussion in Section 7. Here 
we apply the definitions of components given by Tissera et 
al. to our simulations in place of the definitions in Section 3. 

Briefly, in section 2.2 of T12, Tissera et al. define an 
optical radius r op t enclosing 83 per cent of the total galaxy 
stellar mass. A central stellar component is defined by all 
stars with binding energy E less than E ce n , the minimum 
binding energy of stars with r > 0.5r op t. Another threshold 
energy, Fanner, is defined as the minimum binding energy of 
stars with r > 2r op t; stars with E cen < E < Fanner consti¬ 
tute the inner halo component, and stars with E > F?i nner 
the outer halo component 2 . Stars with present-day circular¬ 
ity S > 0.65 are assigned to the disc component, regardless 
of their binding energy. In the inner and outer halo compo¬ 
nents, stars that form within the virial radius of the main 
halo are classified as disc heated if they have S > 0.65 in 
any output of the simulation. Since the disc is ‘hotter’ by 
definition in T12 than in our approach, relatively fewer disc 
stars will be classified as ‘heated’. On the other hand, the 
heated disc category in T12 is not restricted to stars scat¬ 
tered from the z = 0 disc, which could counter the effect 
of a lower circularity threshold by increasing the number of 
stars assigned to this component 3 . 

Stars that are not heated disc are classified as endo de¬ 
bris if they form within the virial radius 4 and debris (which 

2 We find r op t values similar to those of T12. The energy thresh¬ 
olds corresponding to the different components are only slightly 
different from those in fig. 1 of T12: we find (F cen , F?i nner ) = 
(-2.48, -1.56) X 10 5 and (-1.78, -1.11) x 10 5 km 2 S - 2 , for Aq-C 
and Aq-D respectively. If we impose the T12 values of r op t we 
find (F; cen ,F;i nner ) = (-2.22,-1.33) X 10 5 and (-1.90,-1.25) X 
10 5 km 2 s -2 . 

3 Section 2.2 of T13 states that their classification, including the 
identification of the disc, is only applied at redshifts z < 4. The 
C15 results reported in Table Al do not include this restriction. 
If we restrict the assignment of stars to the heated disc category 
to those formed at z < 4 (classifying those formed earlier as endo- 
debris), we find that the only significant change is to the inner 
halo heated disc component of Aq-C, reducing its mass fraction 
from 13 to 6 per cent. 

4 In our implementation of these definitions, membership of the 
self-bound halo identified by SUBFIND is used as a proxy for a 
particle being ‘within the virial radius’. 
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we call accreted) otherwise. Importantly, the endo-debris 
class includes stars that form bound to satellites within the 
virial radius which are later stripped - in the main text we 
classify these as accreted rather than in situ. 

Table A1 compares the masses and mass fractions of 
various components between the two subgrid implementa¬ 
tions and Fig. A1 shows the density profiles of the different 
components. Differences in the properties of the stellar halo 
need to be interpreted in the context of large differences in 
the bulk properties of the central galaxies (summarized in 
Scannapieco et al. 2012). Relative to T12, our galaxies con¬ 
vert about 50 per cent less of their total baryonic mass to 
stars overall 5 but still have more massive discs (by factors of 
6.5 for Aq-C and 1.5 for Aq-D). The central (‘bulge’) com¬ 
ponent is less massive in our simulations (by 49 and 28 per 
cent respectively). The inner halo component differs even 
more - roughly a factor of 10 less massive in our Aq-C and 
a factor of 6 less massive in our Aq-D. The outer halo is 
most closely comparable, being only slightly less massive in 
our simulations. 

Given these differences, it is more useful to compare 
mass fractions than absolute masses. In the Tissera et al. 
inner halo component, we find a substantially smaller contri¬ 
bution from scattered disc stars (discussed in the main text) 
and a much larger contribution from endo-debris. The bal¬ 
ance between in situ and accreted stars is slightly in favour 
of the former in our simulations thanks to this larger endo- 
debris contribution and a much lower mass of accreted stars 
in the same region. As shown in Fig. Al, most of the endo- 
debris in our simulation comes from star formation in galax¬ 
ies after they have become satellites. In the outer halo, the 
masses of endo-debris stars in our simulations are similar 
to those reported by Tissera et al., but the accreted stel¬ 
lar mass in this region is roughly two to three times lower, 
giving a higher endo-debris fraction. Since endo-debris stars 
formed in satellites constitute the majority of endo-debris 
outside ~ 30 kpc in our simulations, it is not surprising that 
the endo debris and debris (accreted) profiles in Fig. Al are 
very similar. 

We conclude that the star formation model in our simu¬ 
lations results in satellites that are more gas rich at accretion 
and that form stars more vigorously thereafter, leading, in 
the language of Tissera et al., to more endo-debris and less 
debris. This is consistent with star formation overall being 
more efficient in the Tissera et al. simulations. 

It is interesting in this regard that our simulations and 
those of Tissera et al. are consistent with the Milky Way 
satellite luminosity function, although this comparison in¬ 
volves only a small number of objects. The gas content of 
surviving satellites and dwarf galaxies outside the virial ra¬ 
dius also provides a useful constraint on predictions for the 
formation of in situ from accreted gas. In our simulations, 
we examine the ratio of neutral (T < 10 4 K) gas mass to 
stellar mass. We find 3, 8 and 9 satellites at r < 250 kpc 
with Mtq/M* > 1 x 10 -4 in haloes C, D and E respectively, 
only slightly in excess of the same count around the Milky 
Way (4) and M31 (5, including M33). Outside this radius, 
we find that our simulated dwarf galaxies have a median 


5 Our simulations would convert slightly more baryons into stars 
at the lower resolution of T12 - see appendix B. 
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Figure Al. Density profies of stellar components in our simula¬ 
tions Aq-C and Aq-D according to the definitions of T12 and T13 
(see the text). Thin lines show different categories of halo stars - 
accreted (dotted; called debris by T13), endo-debris (solid) and 
heated disc (dot-dashed). In addition we show separately the pro¬ 
file of stars that meet the T13 definition of endo-debris but which 
we would classify as accreted, because they form bound to satel¬ 
lites (long dashed). Thick lines show profiles of all stars in different 
phase space ‘regions’ - disc, inner and outer haloes. Arrows mark 
the radius r 0 pt defined by T12. The ‘central’ component of T13 
is confined to log 10 r /kpc < 0.8 and hence is not shown. 


gas fraction of Miq/M* ~ 10, considerably higher than the 
observed median for the Local Group, Mh x /M* ~ 1 (Mc- 
Connachie 2012). This supports the idea that, in our models, 
a low star formation efficiency before infall may lead to high 
star formation rates in stripped gas and in recently accreted 
satellites. 

As noted in the main text, the high star efficiency of the 
central galactic spheroid, presumably fuelled by smoothly 
accreted gas, might also explain the substantially lower 
metallicity of the endo-debris component in T13. Another 
point worthy of note is that our simulations predict even 
lower halo stellar mass from disc heating than Tissera et al., 
despite producing more massive discs at z = 0. 


APPENDIX B: NUMERICAL RESOLUTION 

We have reported three mechanisms responsible for the gen¬ 
eration of in situ halo stars in our simulations - disc heat¬ 
ing, star formation in stripped gas and star formation in 
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Table Al. Columns compare results for haloes Aq-C and Aq- 
D as simulated by Tissera et al. (2012, 2013; T12, T13) and in 
this work (C15), using the definitions given in T12 and T13 in 
all cases. The first row gives the optical radius defined by T12 as 
enclosing a stellar mass M* 5 q pt (second row) equal to 83 per cent 
of the total stellar mass of the galaxy, M*. Rows 3 to 7 give the 
mass of stars in different galactic components of T12, described in 
the text. Rows 8 to 13 separate the inner halo component by origin 
according to T13, giving the mass of ‘debris’, ‘heated disc’ and 
‘endo-debris’ stars and their relative fractions. Rows 14-19 give 
the same breakdown for the outer halo component. All quantities 
assume h = 0.73. 




Aq-C 

T12 

C15 

Aq-D 

T12 

C15 

Topt 

[kpc] 

16.0 

11.0 

14.8 

19.8 

M*, opt 

[10 9 M Q ] 

8.12 

3.49 

6.04 

2.86 

A/disC 

(10 9 M q ] 

2.19 

14.3 

12.2 

18.3 

Accent ral 

(10 9 M q ] 
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32.9 

9.37 

A/inner 

(10 9 M q ] 

19.5 

1.99 

18.6 
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A/outer 

(10 9 M q ] 

4.93 
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5.48 
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Af-total 

(10 9 M q ] 

78.4 

42.0 
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Inner Halo 

A^debris 

(10 9 M q ] 
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0.19 
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0.90 

A^heated 

(10 9 M q ] 
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0.27 
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0.39 
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% 
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9 

43 

31 
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56 
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A/ieated 
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0.07 
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60 
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% 

0 
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0 

2 
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smoothly accreted gas - and discussed similarities and dif¬ 
ferences with other simulations at comparable resolution. All 
three mechanisms are all likely to be sensitive to numerical 
resolution, for different reasons. 

The disc heating dynamics in 7V-body simulations and 
their resolution dependence have been studied extensively 
(e.g. Kazantzidis et al. 2008; Read et al. 2008; Purcell et al. 
2010); in ab initio simulations there is an additional depen¬ 
dence on the subgrid physics governing the formation and 
structure of the disc (e.g DeBuhr et al. 2012; Aumer et al. 
2013). This dependence is now being tackled directly, with 
some success, by recent generations of hydro dynamical sim¬ 
ulations (Scannapieco et al. 2012; Aumer & White 2013; 
Schaye et al. 2015; Snyder et al. 2015). The behaviour of 
subgrid star formation models is less well understood in the 
regime of dwarf galaxies, where there are fewer constraints 
- the most important being comparison to the satellite lu¬ 
minosity functions of the Milky Way and M31 (e.g. Parry 
et al. 2012; Zolotov et al. 2012; Sawala et al. 2013). The 
star formation histories of satellites determine their contri¬ 
bution to the accreted halo and their gas content at accre¬ 
tion, which provides the reservoir for star formation in tidal 
and ram-pressure streams. The formation of stars in gaseous 


streams and clumps depends on the treatment of hydrody¬ 
namic interactions within circumgalactic haloes of hot gas 
(e.g. Mailer & Bullock 2004). Moreover, although these hot 
gas haloes are a robust prediction of the CDM model, their 
detailed properties in simulations are known to be sensitive 
to the treatment of thermal instabilities and galactic winds 
(and, potentially, active galactic nuclei). Subgrid star for¬ 
mation models are typically calibrated in the context of a 
dense cold ISM, with star formation in unstable pockets of 
halo gas most often considered a ‘nuisance’ effect (e.g. Kauf- 
mann et al. 2006; Keres et al. 2012; Parry et al. 2012; Hobbs 
et al. 2013). 

Since our main purpose is to present the phenomenol¬ 
ogy of particular sets of hydrodynamical simulations, rather 
than make predictions for specific observations, we have not 
carried out detailed convergence studies to address these is¬ 
sues. It is possible that the numerical effects given above are 
dominant over physical behaviour in all cases; hence, tests 
of convergence will be absolutely essential for any robust 
predictive model of in situ stellar halo formation. 

In this respect we believe the formation of stars from 
‘smoothly accreted’ hot gas to be the least robust of the 
mechanisms we identify, having three likely sources of reso¬ 
lution dependence: the behaviour of the subgrid star forma¬ 
tion model, the treatment of thermal instabilities and the 
classification of particles as ‘smoothly’ accreted. It is not at 
all clear a priori which of these dependences is most impor¬ 
tant for convergence. 

In Fig. B1 we illustrate the resolution dependence of the 
density profile of stars formed from gas accreted smoothly by 
the main halo. The solid red line corresponds to Aq-C at our 
default resolution level (level 4), including the disc and the 
stellar halo. Solid lines of lighter colours correspond to sim¬ 
ulations from the same initial conditions with progressively 
lower resolution (levels 5, particle mass m p = 2.11 x 10 6 , 
and 6, m p = 1.69 x 10 7 ; these simulations are also discussed 
in Parry et al. 2012). 

As resolution increases (from L6 to L4), the mass of 
stars formed from smoothly accreted gas is reduced at all 
radii. Although the profile of disc stars associated with this 
component appears to converge at resolution L5, this is not 
the case for the stellar halo. In a region r > 30kpc, we find 
4x 10 8 Mq of stars formed from smoothly accreted gas at L6, 
representing 11 per cent of the stellar mass and 0.4 per cent 
of all baryons in the same region. With increasing resolution, 
we see a reduction in the total mass of this component by 
~ 50 per cent for each level, a similar change in its mass 
relative to all stars (from 11 to 4 per cent, a 43 per cent 
reduction in this fraction per level) and a somewhat smaller 
reduction in its mass relative to all baryons (17 per cent 
in this fraction per level). Dashed lines in Fig. B1 tell a 
similar story for the gas that is smoothly accreted but not 
converted to stars - this is always the dominant baryonic 
mass component at r > 30 kpc. 

The absolute change in the mass of halo stars formed 
from smoothly accreted gas is explained only in part by a 
reduction in the total baryonic mass of the halo by ~ 30 
per cent from L6 to L4 (~ 36 per cent for r > 30 kpc; this 
can be compared to a 6 per cent decrease in total mass 
bound to the main substructure). Another significant factor 
must therefore be that the classification of a gas particle as 
‘smoothly accreted’ may also not be converged, because gas 
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Figure Bl. Present-day mass density profiles of smoothly- 
accreted gas particles (dashed) and all star particles formed from 
smoothly accreted gas (solid), in halo Aq-C at our default resolu¬ 
tion level (L4). Thinner lines correspond to lower resolution levels 
L5 and L6 as described in the text. 


associated with DM substructures resolved only at higher 
resolution will be classified as ‘smoothly’ accreted at lower 
resolution. This is illustrated by the results of Wang et al. 
(2011), who examine the total mass (baryons and DM) ac¬ 
creted smoothly on to the main halo in collisionless versions 
of the same initial conditions we simulate. Their fig. 9 shows 
a decrease of ~ 13-15 per cent in the faction of ‘smooth’ ac¬ 
cretion from L5 to L4 (depending on how ‘smooth’ accretion 
is defined). If all DM subhaloes supplied baryons to the main 
halo in the universal ratio, we would expect a similar reduc¬ 
tion in the fraction of gas classified as smoothly accreted. 

We conclude that, at lower resolution, more baryons 
are retained by the main halo and relatively more of these 
baryons are considered as smoothly accreted. The overall 
efficiency with which smoothly accreted gas is converted to 
stars also varies systematically with resolution, although (in 
our simulations) this is a minor contribution to the change 
of in situ stellar halo mass, relative to the first two effects. 



